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A set of Modified Embedded Atom Method (MEAM) potentials for the interactions between Al, 
Si, Mg, Cu, and Fe was developed from a combination of each element's MEAM potential in order to 
study metal alloying. Previously published MEAM parameters of single elements have been improved 
for better agreement to the generalized stacking fault energy (GSFE) curves when compared with 
ab-initio generated GSFE curves. The MEAM parameters for element pairs were constructed based 
on the structural and elastic properties of element pairs in the NaCl reference structure garnered 
from ab-initio calculations, with adjustment to reproduce the ab-initio heat of formation of the most 
stable binary compounds. The new MEAM potentials were validated by comparing the formation 
energies of defects, equilibrium volumes, elastic moduli, and heat of formation for several binary 
compounds with ab-initio simulations and experiments. Single elements in their ground state crystal 
structure were subjected to heating to test the potentials at elevated temperatures. An Al potential 
was modified to avoid formation of an unphysical solid structure at high temperatures. The thermal 
expansion coefficient of a compound with the composition of AA 6061 alloy was evaluated and 
compared with experimental values. MEAM potential tests performed in this work, utilizing the 
universal Atomistic Simulation Environment (ASE), are distributed to facilitate reproducibility of 
the results. 

PACS numbers: 61.50.Lt, 62.20.D- 61. 72. J- 68.35.-p 
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I. INTRODUCTION 

Historically, materials have been developed through 
the correlation of processing and properties. Several im- 
plementations of materials science principles have given 
birth to an engineering framework for materials design. 
Over the past two decades, more efficient computational 
methodologies have been developed and the computa- 
tional power have increased enormously, making the com- 
putational materials design an essential cost-effective tool 
to design materials properties. Since materials complex- 
ities can limit the degree of predictability, several time- 
and length-scale methodologies (hence spatiotemporal 
hierarchy) for computational materials design naturally 
evolved (cf . HorstemeyerS for a review) . Out of several 
computational methodologies, atomistic simulations not 
only can predict the materials properties from a statistical 
viewpoint, but can also quantify the mechanisms of the 
structure-property relationship. One of the most critical 
components of atomistic simulations is the interatomic po- 
tential, which determines the forces on individual atoms. 
First-principles calculations certainly are capable of pro- 
viding very reliable interatomic potentials in a variety of 



chemical environments. However, realistic simulations of 
alloy systems, which are essential to reveal many macro- 
scopic materials properties, often require a number of 
atoms that renders these methods impractical - they ei- 
ther require too much computer memory or take too long 
to be completed in a reasonable amount of time. One al- 
ternative is to use (semi-)empirical interaction potentials 
that can be evaluated efficiently, so that the atomistic 
approaches that use them can, in certain cases, handle 
systems with more than a million atoms. 

The Embedded- Atom Method (EAM) is a widely used 
atomic level semiempirical model for metals, covalent 
materials, and impurities^. MEAM (Modified EAM) in- 
corporates angular dependency of electron density into 
EAM. Atomistic simulations of a wide range of elements 
and alloys have been performed using MEAM potentials. 
MEAM model was first used for silicon, germanium, and 
their alloys^. It was applied to 26 single elemental and to 
silicon-nickeP alloys and interfaces. Gall et aLEihave used 
MEAM to model tensile debonding of an aluminum-silicon 
interface. Lee and Baskes'^ improved MEAM to account 
for the second nearest-neighbor interactions. Also, Huang 
et al^ used MEAM and two other potentials to deter- 



mine defect energetics in beta-SiC. ME AM parameters 
for nickel^ and molybdenum-silicon systeni"'^^' were deter- 
mined by Baskes. MEAM potentials for Cu, Ag, An, 
Ni, Pd, Pt, Al, and Pb based on the first and the second 
nearest-neighbor MEAM were constructed by Lee et alJ^. 
Hu et aL'^2E3 proposed a new analytic modified EAM 
many-body potential and applied it to 17 hep metals. The 
structural properties of various polytypes of carbon were 
described using a MEAM potential by Lee and Leel^. 
Recent work of Lee et al.'^ summarized available MEAM 
potentials for single elements and alloys. Several of these 
potentials were then used to perform large scale atomistic 
simulations to understand the intriguing nature of the 
ductile and brittle fracture^ ^, str ucture-property relation- 
shipi3^ disloc ation dynamica^^E^, and nature of materials 
fractur6pl2i]. 

Aluminum, magnesium, copper, and iron alloys are 
being used in developing materials with novel properties. 
Great popularity of these alloys is connected to their gen- 
eral functional properties, mechanical properties, mass 
density, corrosion resistance, and machinability. Light 
metal alloys, such as magnesium and aluminum alloys, 
are now demanded for use in the automotive and avia- 
tion industries. They performed remarkably well for the 
purpose of decreasing the operating expenses and fuel 
consumption. These alloys usually contain several other 
minor elements, such as silicon, nickel, and manganese, 
and are known to have very complex phase compositions. 
Assessment of such complex systems is a very challenging 
task, since different constituent elements can form differ- 
ent phases, whose selection depends on the ratio between 
the constituents and also on a variety of processing and 
treatment factors. 

Contrary to DFT potentials, most of the single ele- 
ment semiempirical potentials do not combine easily into 
multi-component alloy models. The difficulty of com- 
bining single element EAM potentials into alloy systems 
comes from the need of their normalizatiorr^. The pro- 
cedure to form EAM alloy param eteriz ation from single 
element potentials was suggestecP^^^l^ ^^^ it ^joes not 
guarantee that the resulting potential will be suitable for 
modeling compoundgS^. Alloy potentials usually intro- 
duce new parameters for each pair of elements, allowing to 
fit properties of their binary compounds. The number of 
parameters to adjust and the number of tests to perform 
is proportional to the square of the number of constituent 
elements. In the present MEAM approach, each pair 
interaction is characterized by a total of 13 parameters 
(Table |Vj and the ratio of density scaling factors po for 
constituent elements. Table IT]). Adoption of the default 
value Cniax — 2.8 leads to 9 adjustable parameters for 
each pair. Comparable angularly dependent potentials for 
the Fe-Ni systenPS! also have 9 adjustable pair parameters. 

While the semiempirical potentials have been developed 
and tested for binary alloys^^El]^ binaries, similarly to 
single element potentials, may not combine easily into 
ternaries. Modeling of ternary systems faces a challenge 
since less experimental properties are available for ternary 



systems. Ternary potentials are usually examined only at 
a particular composition range — the number of possible 
compositions grows to the power of the number of con- 
stituent elements. It is also nontrivial to find an equilib- 
rium structure for complex systems of representative size 
at low temperatures. Ternary potentials are only available 
for Fe/Ni-Cr-d22] (MEAM), Pu-Ga-HtPfMEAM), Fe- 
Ti-C/R^Cu-Zr-Ag, Ga-In-N, Fe-Nb-C/I>«Pffl2ll (MEAM) , 
H-C-d2Sl (Reactive Force Field, ReaxFF), Ni-Al-Epl^ Zr- 
Cu-Apl, and Fe-Cu-Np (EAM) systems. To extend from 
binaries to ternaries, MEAM provides a ternary screen- 
ing parameter C'^^^. In the present work we did not 
examine ternary systems. Instead, we performed thermal 
expansion simulations of a compound including all species 
of the potential. The default values of Cmin = 2.0 and 



Cn 



2.8 were applied for ternary screening. Since an 



effort beyond the scope of our project is required for sat- 
isfiable validation of the 5-element alloy potential under 
varying temperatures, compositions, and configuration 
states, we concentrated on basic tests and on providing 
tools to facilitate reproducibility of the testa^. 

In the present study we develop a MEAM potential for 
aluminum, silicon, magnesium, copper, iron, and their 
combinations. We fit the potential to the properties of 
single elements and element pairs, but the model implic- 
itly allows calculations with any combination of elements. 
We show the applied MEAM methodology in Appendix [X 



III 



The DFT calculations are described in Sec. |IIj In Sec 
the single-element volume-energy curves in basic crystal 
structures, and also important material properties, such 
as formation energies of vacancies, self-interstitials, sur- 
faces, and generalized stacking fault energies from MEAM 
are examined and compared with DFT calculations. In 
Sec. |IV[ the MEAM potential parameters for each unlike 
element pair are initialized to fit the ab-initio heat of 
formation, equilibrium volume, and elastic moduli of the 
hypothetical NaCl reference structure. Heat of formation 
of binary compounds in a variety of crystal structures 
from MEAM are thereafter examined and compared with 
the ab-initio and experimental results. The MEAM pa- 
rameters are adjusted to match the DFT formation energy 
of the most stable compounds. The structural and elastic 
properties for several binary compounds and formation 
energies of substitutional defects are compared with ab- 
initio and experimental results. Finally we performed 
thermal expansion simulations of a compound with the 
composition of an AA 6061 alloy (IV C). We conclude 
with a short summary. 



II. AB-INITIO CALCULATIONS 

Ab-initio total energy calculations in this work were 
based on density functional theory (DFT), using the pro- 
jector augmented-wave (PAW) methocpSl as implemented 
in the VASP code^. Exchange-correlation effects were 
treated by the generalized gradient approximation (GGA) 
as parameterized by Perdew et al. ^^ . All DFT calcula- 



TABLE I. Set of the MEAM potential parameters for single elements. The reference structures for Al, Si, Mg, Cu, and Fe 
are fee, diamond, hep, fee, and bee, respectively. Ec is the cohesive energy, ao is the equilibrium lattice parameter, A is the 
scaling factor for the embedding energy, a is the exponential decay factor for the universal energy, l3'-^~^' are the exponential 
decay factors for the atomic densities, f'""^' are the weighting factors for the atomic densities, Cmax and Cmin are screening 
parameters, po is the density scaling factor that is relevant only for element pairs. Definition of these parameters may be found 
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t(0) 
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^(3) 
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^max 


Po 


Al 


3.353 


4.05 1.07 


4.64 2.04 


3.0 


6.0 


1.5 


1.0 


4.50 


-2.30 


8.01 


0.8 


2.8 


1.0 


Si 


4.63 


5.431 1.00 


4.87 4.4 


5.5 


5.5 


5.5 


1.0 


2.05 


4.47 


-1.80 


2.0 


2.8 


2.2 


Mg 


1.51 


3.194 0.8 


5.52 4.0 


3.0 


0.2 


1.2 


1.0 


10.04 


9.49 


-4.3 


0.8 


2.8 


0.63 


Cu 


3.54 


3.62 1.07 


5.11 3.634 


2.2 


6.0 


2.2 


1.0 


4.91 


2.49 


2.95 


0.8 


2.8 


1.1 


Fe 


4.28 


2.851 0.555 


5.027 3.5 


2.0 


1.0 


1.0 


1.0 


-1.6 


12.5 


-1.4 


0.68 


1.9 


1.0 



tions were performed in high precision with the plane- 
wave cut-ofF energy set to 400 eV in order to achieve the 
convergence of heat of formation and elastic properties. 
Integration over the irreducible Brillouin zone was per- 
formed using the F-centered Monkhorst-Pack scheme^ 
with the size gradually increased to 7x7x7 for point 
defects, to 19x19x1 for surfaces, and to 29x29x29 to 
improve convergence of shear moduli at small strains. 
Elastic constants presented here were obtained without 
relaxation of atomic positions. Since most of the examined 
high energy structures are, at best, metastable, relaxation 
does not maintain the crystal symmetry, resulting in large 
energy changes and unphysical elastic constants. 



of DFT energy differences and volume ratios can possi- 
bly be obtained by optimization of Si and Cu MEAM 
parameters. Fe MEAM potential applied in the present 
work is a MEAM-p variant of Fe potential from the recent 
effort of Lee et al.'^, exhibiting a correct low tempera- 
ture phase stability with respect to the pressure. The fee 
equilibrium energy and volume from this Fe potential is 
very close to the bcc equilibrium in order for the struc- 
tural transition to appear at finite temperature without 
magnetic contribution. In general the MEAM potentials 
of the present work reproduced the DFT results for the 
individual elements fairly well. 



III. MEAM PARAMETERS FOR SINGLE 

ELEMENTS 

The present MEAM parameters for single elements are 
listed in Table |I] The initial values of the se param eters 
were taken from existing MEAM potentialJ^ElETlill. The 
Cniin screening parameter for Al, Mg, and Cu was lowered 
from 2.0 to 0.8 to improve the GSFE curves (Sec. |IIIE |. 
The Mg potential was adjusted to reproduce the DFT 
values of hep, bcc, and fee energy differences, vacancy 
formation energy, and (lOlO) surface formation energy. 
The Al potential was modified to prevent formation of an 



unknown structure at elevated temperatures (Sec. IV C I 



A. Energy dependence on volume of single 

elements in fee, hep, bcc, and simple cubic crystal 

structures 

The first test of the validity of MEAM potential for sin- 
gle elements is a comparison of the energy-volume curves 
in the fee, hep, bcc, diamond, and simple cubic crystal 
structures, shown in Fig. [T] The MEAM potentials ap- 
propriately capture the lowest energy structures of Al 
(fee). Si (dia), Mg (hep), Cu (fee), and Fe (bcc). Also, 
the equilibrium volumes of several crystal structures from 
MEAM closely match the DFT results. Better match 



B. Vacancies 

The formation energy of a single vacancy EJ^'^ is defined 
as the energy cost to create a vacancy: 

Er = Etot[N]-Ne, (1) 

where Etot[N] is the total relaxed energy of a system 
with TV atoms containing a vacancy and e is the energy 
per atom in the bulk. Cell volume and atomic positions 
were relaxed in each case. Table [III shows the formation 
energies of a single vacancy for the fee Al cell, diamond 
Si cell, hep Mg, fee Cu, and bcc Fe obtained from the 
MEAM and DFT calculations. The MEAM systems sizes 
were 5x5x5 primitive fee and bcc cells, 3x3x3 primitive 
diamond cells, and 8x4x4 orthogonal hep cells. For 
the DFT systems, the simulation sizes were 5x5x5 fee, 
4x4x4 diamond and bcc, and 4x4x2 hep primitive cells. 

The vacancy formation energy of Mg was slightly im- 
proved in comparison with previous MEAM results^. 
Overall agreement of vacancy formation energies between 
MEAM, experiment, and DFT was within a few eV, and 
the present results are comparable or better than those 
from other calculations. The reduction in volume due 
to the formation of a vacancy agrees well with the DFT, 
except the value for Fe is somewhat low. 
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FIG. 1. (Color online) Energy-volume dependenee of Al, Si, 
Mg, Cu, and Fe in fee, hep, bee, diamond, and simple cubic 
crystal structures relative to the ground state. 



TABLE II. Calculated single vacancy properties. Single va- 
cancy formation energy EJ^'^ and formation volume Q^ values 
are obtained from the relaxed structures containing single va- 
cancies. Here f2o is the bulk atomic volume. All energy values 
are listed in eV. The results from the MEAM calculations are 
compared with the results from the DFT calculations given 
inside the parentheses, other simulations, and experiments. 



a„/no 



Al 
Si 

Mg 
Cu 
Fc 



Present 
0.67 (0.5) 
3.27 (3.6) 
0.89 (0.7) 
1.10 (1.0) 
1.65 (2.1) 



_Othcrs 
I.6^0.5J 

k6H3.q 




Present 

0.67 (0.7) 

0.21 (0.3) 

0.72 (0.8) 

0.75 ("-0.9) 

0.47 ('-O.S) 




"^ MEAM results by Lee et al^ 

'' Calculated using EAM parameters extracted from Liu et al)^ 
" DFT calculation by WriglitEHl 
d AMEAM results by Hu et ai.HS 
° Experimental results by Tzanetakis et aZ.EIl 
' DFT results by Carling et al^ 
s Experimental value by HehenkampEU 
^ DFT calculation by Andersson and SimakE31 
' Experimental value by Dannefaer et aLl^, Throwe et a/.l^ 
■i Experimental value by Hehenkamp et a^ESl 
^ DFT value by Krimmel and FahnleES 
' EAM value by Mendelev et alW^ 
™ Finnis-Sinclair potential value by Ackland et a^ESl 
" DFT value by Domain and Becquartl^ 
° Experimental value by Schaefer et aLED 
P Experimental value referenced in Ackland et ai.ESl 
^ EAM value by Mendelev et a/.EH 
"■ EAM value by Mishin et alEU 
= MEAM value by Timonova^ aZ.Ell 
* MEAM value by Ryu et al. ^^ 



C. Self-interstitials 

The formation energy of an interstitial point defect E"™* 
is given by 



Ej"* 



Etot[N + l]~Nex-SY 



(2) 



where i^tot [N + 1] is the total energy of a system with 
N type-X bulk atoms plus one impurity atom of type-Y 
inserted at one of the interstitial sites, and £x (£y) is 
the total energy per atom of type-X (type-Y) in its most 
stable bulk structure. The inserted atom Y can be of the 
same type as the bulk, in which case the point defect is 
called a self- interstitial defect. Self-interstitial formation 
energies were calculated for Al, Si, Mg, and Cu at the 
octahedral, tetrahedral, and dumbbell sites. Dumbbell 
orientations were [100] for fee, [0001] for hep, and [110] for 
bcc and diamond structures. Relaxations of the atomic 
positions and the volume were also performed, and the 
DFT and MEAM results are listed in Table UHl Similar 
to the previous calculations, the MEAM systems sizes 
were 5x5x5 primitive fee and bcc cells, 3x3x3 primitive 
diamond cells, and 8x4x4 orthogonal hep cells. For 
the DFT systems, the examined sizes were 5x5x5 fee 
primitive cells, 4x4x4 diamond and bcc primitive cells, 
and 4x4x2 hep primitive cells. 



TABLE III. The formation energies of various Al, Si, Mg, Cu, 
and Fe self-interstitials. All energy values are given in eV. The 
results from the MEAM calculations are compared with the 
DFT results and other classical MD (CMD) simulations. 



Interstitial 


DFT MEAM 






CMD 




DFT 


Al 


Present 


RefP RefP 


Ref.E3 


Ref.l^ 


tetrahedral 


3.3 3.32 


3.1f 


a 


2.94 






octahedral 


2.8 3.26 


3.0f 


a 


2.82 






split (100) 


2.7 2.77 


2.6J 


a 


2.46 


2.59 


2.43 


Si 


Present 


Ref.«« 


Ref.^^ 


Ref.Ea 


Ref.E2l 


split (110) 


3.7 3.71 


3.88 


4.7 


3.9 


3.40 


Mg 


Present 


Ref.E2l Ref.l32l 


Rei^ 


Ref.E3 


tetrahedral 


2.2 1.63 


1.5^ 


a 


1.53 


1.41 


2.35 


octahedral 


2.2 1.57 


2. If 


a 


1.29 


1.20 


2.36 


split (0001) 


2.3 1.78 


1.5^ 


a 








Cu 


Present 


Ref!^ 


Reim 


Reim 


Reim 


tetrahedral 


3.9 3.37 


2.9; 


a 








octahedral 


3.5 2.72 


2.9'; 


a 








split (100) 


3.3 2.59 


2.81 


3.06 


3.23 


2.93 


Fe 


Present 


Ref.l^ 


1 Ref.ESl 


Ref.^ 


Ref.^ 


tetrahedral 


4.2 4.31 


4.1f 


a 






4.14 


octahedral 


5.0 4.78 


4.1C 


a 






4.82 


split (110) 


3.9 3.79 


3.53 


4.11 


3.65 


3.64 


split (111) 


4.9 4.28 


4.02 


4.01 


4.24 


4.34 


split (100) 


4.8 4.81 


4.3 


1 


4.28 


4.60 


4.64 



Calculated using parameters from the given reference. 



In general, the DFT results are well reproduced or 
slightly underestimated by the MEAM potentials. Ac- 
cording to the present MEAM potential, the most stable 
form of a self-interstitial defect for fee Al is a dumbbell 
along the [100] direction, in agreement with the DFT re- 
sults and an experimental observation by Jesson et alM^. 
The results for Mg are better than those published previ- 
ously"*^ ''^ . The present Mg potential indicates that the 
tetrahedral site will be most stable in agreement with 
the DFT calculations. For both Cu and Fe, the new 
MEAM potential produces the same relative stability of 
the examined interstitial sites with the DFT calculations. 



TABLE IV. Surface formation energies for Al, Si, Mg, Cu, and 
Fe. The units are mj/m^. The second column indicates if 
the structure was relaxed. Comparisons with other classical 
MD (CMD), DFT, and experimental values for polycrystalline 
surfaces and Si facets are also given. 



Surface Rlx DFT MEAM 



CMD 



DFT EXP 



Al 




Present 


Ref.E] 


Rei^^ 


Ref.E^ 


Ref.ES 


(111) 


No 


780 


820 




913 


1199 


1143 


(111) 


Yes 


780 


752 


629 


912 






(110) 


No 


990 


1154 




1113 


1271 




(110) 


Yes 


960 


1135 


948 


1107 






(100) 


No 


890 


1121 




1012 


1347 




(100) 


Yes 


890 


1088 


848 


1002 






Si 




Present 


Ref.l^ 


Ref.l^ 


Ref.l^ 


Ref.l^ 


(111) 


No 


1620 


1254 


1405 




1820 




(111) 


Yes 


1570 


1196 


1405 




1740 


1230 


(100) 


No 


2140 


1850 


2434 




2390 




(100) 


Yes 


2140 


1743 


1489 




2390 




(100) 


2x1 




1241 




2050 


1450 


1360 


Mg 


Rlx 


Present 


Ref.EIl 


ReiW^ 


Ref.E^ Ref.E^ 


(0001 


) No 


530 


780 




500 


792 


785 


(0001 


) Yes 


530 


713 


310 


499 






(1010 


) No 


850 


878 




629 


782 




(1010 


) Yes 


850 


859 


316 


618 






Cu 




Present 


Rei^ 


Rei^f' 


Ref.l^ 


Ref.l^ 


(111) 


No 


1290 


1411 


1185 


919 


1952 


1825 


(111) 


Yes 


1290 


1411 


1181 


903 






(110) 


No 


1550 


1645 


1427 


1177 


2237 




(110) 


Yes 


1510 


1614 


1412 


1153 






(100) 


No 


1440 


1654 


1291 


1097 


2166 




(100) 


Yes 


1430 


1653 


1288 


1083 
Ref.l^ 


Ref.I^ 




Fe 




Present 


Ref.^ 


Ref.l^ 


(111) 


No 


2760 


1366 


1941 


2012 


2660 


2475 


(111) 


Yes 


2700 


1306 


1863 


1998 


2580 




(110) 


No 


2420 


1378 


1434 


1651 


2380 




(110) 


Yes 


2420 


1372 


1429 


1651 


2370 




(100) 


No 


2500 


1233 


1703 


1790 


2480 




(100) 


Yes 


2480 


1222 


1690 


1785 


2470 





^ Value from the given reference. 

^ Calculated using parameters from the given reference. 



D. Surfaces 

A semi-infinite surface is one of the simplest forms of 
defects. To test the transferability of the new MEAM 
potentials, formation energies for several 1x1x7 surface 
slabs with 8 A vacuum layer were computed. Eight atomic 
layers were used for the Si(lll) surface and 12 layers 
for the 2x1 Si(lOO) surface reconstruction. The surface 
formation energy per unit surface area E'surf is defined as 



T^surf 



Etot[N]-Ne 
A 



(3) 



where i^tot [N] is the total energy of the structure with 
two surfaces, N is the number of atoms in the structure, 
e is the total energy per atom in the bulk, and A is the 
total area of both surfaces. Table IIVI shows the surface 
formation energies of several surfaces constructed from 
fee Al, hep Mg, fee Cu, and bcc Fe crystals. Results 
from the present MEAM potentials agree, in the order 
of magnitude, with the DFT calculations, except for Fe 
values being underestimated. 

The surfaces with lowest energy without reconstruction 
are identified correctly by the present MEAM potentials. 
The 2x1 reconstruction of the Si(lOO) surface leads to 



symmetric dim.ers in accord with other Si potentials^. 
Note that surface formation energies from the present 
PAW GGA calculations are lower than our previously 
published results using ultrasoft pseudopotentials within 
local density (LDA) approximation — it is known that 
GGA leads to surface energies which are 7-16% lower 
than LDA values fo r jelliu m and 16-29% lower than the 
experimen tal re sult d^^^ l ^^ l A procedure^ and new DFT 
functionala^^ESl were suggested to correct the errors of 
LDA and GGA approximations. Similar correction can 
be applied to vacancy formation energiea^Sl, but such 
corrections were not applied in the present study. 



E. Stacking faults 

Using an assumptio n of a planar dislocation core, the 
Peierls-Nabarro modepIEIl jg a powerful theory to quan- 
tify the dislocation core properties. In that model, a dis- 
location is defined by a continuous distribution of shear 
along the glide plane, and the restoring force acting be- 
tween atoms on either sides of the interface is balanced 
by the resultant stress of the distribution. As shown 
in the recent study of Carrez et alJ^, a solution of the 
Peierls-Nabarro model can be obtained numerically by 
identifying the restoring force to the gradient of the gener- 
alized stacking fault energy (GSFE) curve^**. In addition. 
Van Swygenhoven et al.^ claimed that the nature of slip 
in nanocrystalline metals cannot be described in terms 
of an absolute value of the stacking fault energy, and a 
correct interpretation requires the GSFE curve, which 
shows the change in energy per unit area of the crystal 
as a function of the displacement varied on the slip plane. 
However, the GSFE curve is not experimentally accessible. 
Therefore, to model dislocation properties reliably, the 
GSFE curve calculated with the MEAM potential must 
reproduce the DFT data. 

The stacking fault energy per unit area of a stacking 
fault Ef is defined as 



^f ^ A ' 



(4) 



where -Etot [N] is the total energy of the structure with a 
stacking fault, N is the number of atoms in the structure, 
e is the total energy per atom in the bulk, and A is the 
total area of surface. 

As a validation test of the MEAM potential, the GSFE 
curves obtained by molecular statics (MS) were compared 
with the DFT data by Zimmerman et a/.^SI for Al and Cu, 
by the present authors for Fe, and by Datta et al.^ for 
Mg. After lowering the Cmin parameter to 0.8, the GSFE 
curves calculated by MS using the MEAM potential for 
Al, Cu, and Mg show the skewed sinusoidal shape in 
agreement with the DFT predictions (Fig. [2]) illustrating 
reasonable agreement with the DFT GSFE curves. 



IV. MEAM PARAMETERS FOR ELEMENT 
PAIRS 

The MEAM potential parameters for each element pair 
were initialized to match the ab-initio heat of formation, 
equilibrium volume, bulk modulus, and elastic moduli 
in the hypothetical NaCl reference structure, which was 
chosen for its simplicity. Since the equilibrium volume, 
cohesive energy, and bulk modulus of the NaCl struc- 
ture are directly related to MEAM parameters, they can 
be reproduced exactly. An improved agreement of the 
shear moduli from MEAM and ab-initio simulations was 
achieved in some cases by adjusting the electron density 
scaling factor pQ. Then, heat of formation of binary com- 
pounds in a variety of crystal structures from MEAM 
were examined and compared with the ab-initio results. 
To correlate the MEAM results with the lowest forma- 
tion energies of the compounds from DFT calculations, 
the MEAM screening and AiJg^ parameters for element 
pairs were adjusted. The final MEAM parameters are 
given in Table [Vj The predicted MEAM properties for 
the NaCl reference structure are compared with DFT 
results in Table \VT\ and show that in general the MEAM 
heat of formation, bulk modulus, and equilibrium volume 
reproduce the DFT results well. In contrast, the shear 
elastic constants are not well reproduced. In fact the sign 
of the shear elastic constant, representing crystal stability, 
is frequently in disagreement with the DFT results. This 
is really not a significant problem as the NaCl structure 
does not exist in nature. A more important criteria for 
success of these potentials is how they perform for lower 
energy crystal structures. We address this issue in the 
next section. 



A. Heat of formation for binary compounds 

The alloy phases that the MEAM potential predicts as 
most likely to form at the temperature T = K are those 
with the lowest heat of formation per atom, AH, which 
is defined as 



AiJ = 



EtotWx + Ny] - Nxex - Nyey 
Nx + Ny 



(5) 



Etot is the total energy of the simulation cell, TVx and Ny 
are the numbers of type-X and type-Y atoms in the cell, 
£x and £y are the total energies per atom for type-X and 
type-Y in their ground state bulk structures, respectively. 
To check the validity of our new potentials, we com- 
puted the heat of formation per atom for many inter- 
metallic phases of all alloy pairs. The total energy values 
in Eq. ^ for Bl, B2, B3, CI, C15, DO3, A15, LI2, and 
other relevant structures were evaluated at the optimal 
atomic volume for each structure. Heat of formation for 
basic binary compounds based on the new MEAM poten- 
tial and DFT results were calculated and compared with 
experimental values (Figures [sHs]) . The DFT and MEAM 
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FIG. 2. GSFE curves for Al, Mg, Cu, and Fe obtained with the MEAM potential and compared with the DFT data. 

TABLE V. The MEAM potential parameters for element pairs. AH^Y is the heat of formation of the NaCl structure (reference) 
with the type-X and type-Y elements relative to the energies of elemental X and Y in their equilibrium reference state, r^ 
is their equilibrium nearest neighbor distance, a is the exponential decay factor for the universal energy, Cmax and Cmin are 
scree ning pa rameters {C^^'^ denotes type-Y element between two type-X elements). Non-zero parameters Sr = Sa = O-l in Rose 
Eq. (Bl B4| were used for SiFe pair. 
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results for the phases with lowest AH are also shown in 

Tables EmUIXl 

The agreement between MEAM and DFT is quite sat- 
isfactory. In most cases, the MEAM results preserve the 
order of stability predicted by the DFT results. The 
differences in the heat of formation per atom from the 
MEAM and DFT results are less than 0.5 eV at most. In 
general the atomic volumes predicted by MEAM agree at 
least qualitatively with the DFT and experimental results. 
The MEAM calculations of the bulk moduli also agree 
semi-quantitatively with DFT and experimental results, 
usually within 20%. Predicted shear moduli usually fol- 
low the DFT and experimental results, but in some cases 
there is significant disagreement. 



where Etot[{N — 1) -I- 1] is the total energy of a system 
of iV — 1 host type-X atoms and one type-Y atom that 
replaced type-X atom in the original bulk position, ex and 
£y are the total energies per atom for type-X and type-Y 



atoms in their ground state bulk structures. Table VII 



shows the results of substitutional defect calculations us- 
ing the MEAM potentials and the DFT results. In general 
the MEAM results qualitatively agree with the DFT re- 
sults. In a number of cases of small heat of formation, 
MEAM indicates a small heat, but of the incorrect sign. 
The most significant error is for Al in Si where MEAM 
predicts a large endothermic heat and DFT predicts a 
much smaller value, otherwise there is general agreement. 



B. Substitutions 



C. Finite temperature tests 



The formation energy of a substitutional point defect 
in the case of the substitution of a type-X atom of 
the host with a type-Y atom, is defined by 



TT^sub 



psub 



£;tot[(^-l) + l]-(iV-l)ex-eY (6) 



Real life applications of MD potentials require extensive 
testing at finite temperatures. Basic finite temperature 
tests of the potentials, in accord with recommendations 
of Lee et al.'^, revealed formation of an unknown solid 
structure when the temperature of fee Al crystal was 
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FIG. 3. (Color online) Heat of for mati on of AlSLAlMg, AlCu, AlFe binary compounds from MEAM, DFT, and experiments. 
References: AISP, AIM^S, AlCrPSSU, A1F«P^ES11. DFT poin ts are labeled on the left, MEAM and experimental on the right. 
Values for the most stable compounds are also shown in Table IVlll] The inside plot is a magnified portion of a larger plot. 
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TABLE VI. Structural and elastic properties of element pairs 
in the reference NaCl (Bl) crystal structure from DFT and 
MEAM calculations. AH is the heat of formation in eV/atom, 
Vo is the volume per atom in A''. Elastic constants Bq, C44, 
and (Cii — Ci2)/2 are in GPa. 





pair 


method 


AH 


Vo 


Bo 


C44 


C11-C12 
2 


AlSi 


DFT 
MEAM 


0.28 
0.28 


17.9 
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10 
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76 
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AlMg 


DFT 
MEAM 


0.42 
0.23 


23.7 
23.6 


30.9 
33.9 


-18 
-3 
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SiMg 
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SiFc 
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DFT 
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DFT 
MEAM 



DFT 
MEAM 

DFT 
MEAM 
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0.36 
0.26 



16.1 
16.2 

14.7 
14.7 
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77.4 
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92.7 



0.41 
0.20 
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54.9 



0.39 
0.14 
0.25 
-0.07 



14.9 
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MgFe 


DFT 
MEAM 
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17.7 
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50.4 
56.5 
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CuFe 


DFT 
MEAM 


0.78 
0.63 


14.1 
14.2 


107.4 
111.8 


-23 
10 
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131 



TABLE VII. The formation energies of substitutional point 
defects in Al, Si, Mg, Cu, and Fe. All energy values are given 
in eV. DFT values are given in parentheses. 



Host - 



Substitute atom 



Al 



Si 



Mg 



Cu 



Fe 



Al 



0.5 (0.5) -0.2 (0.05) -1.1 (-0.1) -1.3 (-0.4) 



Si 


7.0 (0.9) 
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1.9 (2.) 
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(1.9) 


Mg 
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Cu 


0.7 (-0.7) 
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Fe 


0.2 (-0.7) 


-2.9 (-1.1) 


0.8 (1.0) 
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increased to 800 K under zero pressure conditions. To 
prevent formation of this structure, /3^^' and t^^' parame- 
ters of Al were adjusted. Heating of other elements under 
zero pressure conditions did not result in forming new 
structures. 

To test a system including all components of the new 
potential, an 20-100 °C average thermal expansion coeffi- 
cient of a model system with the composition similar to 
AA 6061 alloy ( Table pC| was evaluated and compared with 
experimental data. Atoms of constituents were placed in 
the substitutional positions of a 20x20x20 fee Al cell. The 
system was heated from -200 °C to 20 °C (and 100 °C) 
over the interval of 0.1 ns, and then equilibrated at 20 °C 



11 



TABLE VIII. Structural and elastic properties of element pairs in varying crystal structures from the present DFT and MEAM 
calculations compared with references and measured values. Aiif is the heat of formation in meV/atom, Vb is the volume per 

Ci2)/2 are in GPa. 
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[Toaf 



Jii5iPSr 

031 



( J125I 



AlFc 



B2 



DFT 



MEAM/EAM 



EXP 



12 



ll.„,j^...„..^--^^ 
11.3^3 11.2^3 
11.6JIq 11.8SS 
ll.g^iMl, 12.1 



174. 
183. 
155.' 
156. 



138. 
107 



61. 
38. J 



-342. -29JS^. -26(P3, 
-10CP3 



12 
13.9 



12.3 



,(P3. -42ii321, _25(P3, 



2isa 



jT^ 



11.41 



145 
138 
193. 
152. 
136. 



.„ 144.(PS1, 

T041 



111.7. 117. d^^ 79.6. 18. (P^. 
110.^ 7.cjml 



1341 

T35I' 



127.1 



43 



(pS, 33.'P3 



DO3 DFT -202. -20.'4^. -20rt^ 11.79. 11 

-221^251. -23fP3. -2nJI2Il. 14.6E@ 

_20llSll -l86Cpa 12.0cIiSI 

MEAM/EAM 346. -20(P1, _22Ji23 12.01. 11 

-7423 12.8ffi21, 

EXP -2oinn _32i[iisi i2.o'in33i 



_J74.C>ii2I, 

1241 



127. 1^^^ 

140.0, 137.Si23J 25.5, lO.di^IT 



A1FC3 



AlFe2 
AI3F0 



AI2F0 



,8llii21_i2.oiLi2(a 

11.8^ 
11.5'PSr 

,7tE21 

10.81 




129. 0. 16 2. (P^ 
126.$^ 

isi.'iEasi 



30.nM.Cp3, 

12.P3 

20.JI3I1 



L12 


DFT 


-1^-ts^^.od^, 


12.13,12.l4ii2J 
U.l&i, 12.3lP231 




166.3. 158. (P^^. 
168.(P^ 


125.1 


11.5 




MEAM/EAM 


274, -18lii22I 


12.88 






139.5 


96.9 


84.0 


A15 


DFT 


-161 


12.08 






156.9 


67.3 


135.1 




MEAM/EAM 


205 


12.59 






166.7 


35.5 


153.5 


C15 


DFT 


115. 6di2il 


12.42, 


11.4ii2iJ 




130.2 


55.0 


52.2 




MEAM/EAM 


925 


13.67 






127.8 


259.2 


551.0 


A15 


DFT 


-161 


13.91 






121.5 


67.7 


120.0 




MEAM/EAM 


321 


15.03 






103.5 


1.8 


66.6 


LI2 


DFT 


122, 15(*i2iJ, lO^iSU 


13.68, 


13.0Tli2iJ, I3.6<>i211 


126.5, 98.di2il 


85.9 


49.9 




MEAM/EAM 


-49 


14.83 






108.5 


59.4 


20.3 


DO3 


DFT 


_25, -9cP^, -ldi241 


13.38, 


13.5#^, 13.3ii211 


126.0, 119.di2il 


91.4 


-48.3 




MEAM/EAM 


266 


14.81 






93.8 


57.7 


-31.1 


Cllb 


DFT 


371, 42di2iJ 


12.78, 


12.3^i2il, 12.8(ii^ 


149.0 








MEAM/EAM 


106 


14.71 













CI 



DFT 

MEAM/EAM 



-72 
-76 



15.25 
16.12 



98.6 
90.4 



76.8 
47.5 



55.0 
36.4 
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TABLE IX. Structural and elastic properties of element pairs in varying crystal structures from the present 
calculations compared with references and measured values. Aiif is the heat of formation in meV/atom, Vb 

Ci2)/2 are in GPa. 



atom in A"^, and elastic constants Bo, C44, and (Cn 



DFT and MEAM 
is the volume per 



compos. 


str. 


met. 


AH 






Vo 


Bo 


C44 


C11-C12 
2 


SiMg2 


CI 


DFT 

MEAM/EAM 

EXP 


^185, -18C^i2U 
42 

-22gnzi 




21.41 
23.05 


54.1 
47.8 


47.6 
20.9 


47.3 
32.2 


SiMg3 


LI2 


DFT 
MEAM/EAM 


-11 

24 






19.29 
20.70 


50.8 
57.5 


29.9 
23.5 


37.1 

21.8 




A15 


DFT 
MEAM/EAM 


69 

-14 






20.09 
21.25 


44.1 
56.3 


9.3 
32.2 


31.8 
33.9 


SiCu3 


LI2 


DFT 

MEAM/EAM 

EXP 


-22, 35 
28 


106 




12.18 
13.37 


137.3 
134.9 


65.0 
64.4 


38.1 
33.0 


SiCu2 


CI 


DFT 
MEAM/EAM 


60 

-41 






14.26 
14.81 


111.9 
102.8 


76.3 
97.9 


23.1 
16.2 


SiFc 


B20 


DFT 

MEAM/EAM 

EXP 


-489, -4842iJ, - 
-132 
4l(P3 


-42C(i22l 


11.04 
13.11 


226.5 








B2 


DFT 
MEAM/EAM 


-457 
-222 






10.55 
13.09 


231.9 

177.7 


87.0 
36.2 


155.8 
225.3 


SiFc3 


DO3 


DFT 

MEAM/EAM 

EXP 


-305, -31^, - 
-269 


28CP2ZJ 


10.99 
12.03 


204.5 
169.2 


142.4 
91.6 


54.5 
36.6 




A15 


DFT 

MEAM/EAM 


-251 
-232 






11.44 
12.28 


173.3 
190.1 


72.3 
47.6 


125.3 
119.9 




LI2 


DFT 
MEAM/EAM 


-236, - 
-149 


.23daz) 




11.38 
11.80 


188.4 
188.6 


116.0 
65.1 


26.6 
135.8 




DO19 


DFT 
MEAM/EAM 


-232, - 
-216 


-23d2Il 




11.28 
12.05 


160.1 






SiaFc 


CI 


DFT 
MEAM/EAM 


-249, - 
-140 


.24dazl 




13.05 
15.53 


169.4 
158.7 


136.8 
95.7 


15.5 
41.3 




C16 


DFT 
MEAM/EAM 


-251, - 
111 


.24*iZl 




12.39 
13.64 


170.6 






SiFe2 


CI 


DFT 
MEAM/EAM 


-12 

-516 






12.95 
12.95 


159.6 
182.4 


82.6 
154.6 


-1.4 
226.8 


MgaCu 


Ci, 


DFT 
EXP 


-129. -13^11. - 
-9cJi^ 


13ii22J, 


18.13 


57.1 







MgCu2 C15 DFT 



MEAM/EAM 
EXP 



-157. -I6d21l. -16d^flS]. 14.58 

-140 14.50 

-ll'itilll, -s^m 



90.6 
104.7 


45.4 
21.5 


25.2 
-12.0 


69.3 
73.6 


60.3 

51.1 


16.9 

-11.6 


96.6 
103.9 


62.6 

42.5 


21.5 
29.0 



MgCu 
MgCu3 



B2 

LI2 



DFT 

MEAM/EAM 
DFT 
MEAM/EAM 



-117 
-28 
-71 
-25 



15.65 
15.56 
13.49 
13.81 



DO3 



DFT 
MEAM/EAM 



-54 
107 



13.52 
13.43 



96.5 
95.7 



75.2 
64.1 



0.3 

-22.4 



MgFc2 


C15 


DFT 
MEAM/EAM 


62 
471 




MgFc3 


LI2 


DFT 
MEAM/EAM 


181 
243 




Mg3Fc 


LI2 


DFT 
MEAM/EAM 


204 
249 




MgFc 


B2 


DFT 
MEAM/EAM 


347, 
238 


35^ 


CuFea 


LI2 


DFT 
MEAM/EAM 


133 

267 





13.58 
14.79 



91.4 
95.8 



72.0 
322.9 



53.4 
605.7 



13.11 
13.65 
18.29 
19.08 



122.6 


96.9 


14.0 


116.4 


64.5 


45.5 


51.6 


52.1 


18.9 


54.5 


30.0 


23.1 


71.4 


68.7 


-25.4 


86.4 


63.2 


13.7 


132.5 


99.9 


6.8 


172.6 


80.5 


68.0 



15.39 
15.82 
11.73 
11.78 





A15 


DFT 


178 


11.95 


137.6 


54.2 


134.7 






MEAM/EAM 


129 


12.26 


192.3 


70.5 


155.2 


CU3FC 


DO3 


DFT 


175 


12.10 


139.2 


105.5 


-2.0 






MEAM/EAM 


287 


11.59 


134.0 


130.8 


20.7 




LI2 


DFT 


224, 34d2ZI 


12.16 


135.2 


60.7 


-0.0 






MEAM/EAM 


230 


12.01 


153.2 


66.7 


71.7 


CuFo 


B2 


DFT 


264, 26^211 


11.88 


211.4 


108.7 


-52.0 






MEAM/EAM 


217 


12.10 


161.7 


105.6 


45.7 
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TABLE X. Composition limits of AA 6061 allojISHl and a 
model system used to estimate thermal expansion coefficient. 







Limits 






Model [wt. 


%] 


Element 


low [wt. ^ 


^o] high [wt. 


%] 


Si 


0.40 




0.8 




0.51 




Mg 


0.8 




1.2 




1.00 




Cu 


0.15 




0.40 




0.30 




Fe 


no 




0.7 




0.50 




Mn 


no 




0.15 




0.00 




Cr 


0.04 




0.35 




0.00 




Zn 


no 




0.25 




0.00 




Ti 


no 




0.15 




0.00 





TABLE XL Thermal expansion coefficient of single crystal Al 
and AA 6061 alloy between 20°C and 100°C in the units of 
/im/m/K. 





CMD 


Exp 




present Ref.li^ 


Ref.U421 RgflliOl 


Alfcc 
A A 6061 


14.4 15-25 
14.6 


25.4 23.6 
23.6 



(and 100 °C) for 1 ns under zero pressure conditions. 
Table pCl| shows the values of 20-100 °C average ther- 
mal expansion coefficients. The MEAM result for single 
crystal Al is in the lower range of other MD potentials 
and experiments. Since Al is a dominant element of 
the AA 6061 alloy, the thermal expansion coefficient for 
alloy is similarly underestimated, possibly also due to 
imperfections of the structure of the real material. 



V. CONCLUSIONS 

In this study we developed MEAM potentials for the 
pair combinations of aluminum, silicon, magnesium, cop- 
per, and iron. The MEAM formalism allows any of these 
potentials to be combined to enable prediction of multi- 
component alloy properties. These potentials reproduce a 
large body of elemental and binary properties from DFT 
calculations at the temperature of K and experimental 
results. Basic finite temperature tests of the single ele- 
ment potentials and their alloy combinations were also 
performed. With focus to facilitate reproducibility of 
the presented results^, and subject to further testing 
and improvements, these potentials are one step towards 
designing multi-component alloys by simulations. 



porting this study. Computer time allocation has been 
provided by the High Performance Computing CoUabora- 
tory (HPC^) at Missi ssip pi State Univ ersity. Computa- 
tional package LAMMP&^^with ASlP^ interface was used 
to perform MD simulations. Much appreciated tests of the 
new MEAM potentials, including the high temperature 
simulations of Al that revealed formation of unknown Al 
phase at 800 K, were performed by Chandler Becker and 
Tanner Hamann at the Metallurgy Division of the Material 
Measurement Laboratory, National Institute of Standards 
and Technology (NIST) . Comparison of ab-initio elastic 
constants and related discussion with Hannes Schweiger 
from Materials Design are also appreciated. Classical MD 
potentials from other authors examined in this study were 
downloaded fro m th e Interatomic Potentials Repository 
Project databaseHni. 

Sandia National Laboratories is a multi-program lab- 
oratory managed and operated by Sandia Corporation, 
a wholly owned subsidiary of Lockheed Martin Corpora- 
tion, for the U.S. Department of Energy's National Nu- 
clear Security Administration under contract DE-AC04- 
94AL85000. 



Appendix A: MEAM theory 

The total energy E oi& system of atoms in the MEAMS 
is approximated as the sum of the atomic energies 



E = Y^E,. 



(Al) 



The energy of atom i consists of the embedding energy 
and the pair potential terms: 



E, - F, {p. 






4>tj {rtj) ■ 



(A2) 



F is the embedding function, pi is the background electron 
density at the site of atom i, and (pij (rij) is the pair 
potential between atoms i and j separated by a distance 
rij. The embedding energy Fi {pi) represents the energy 
cost to insert atom i at a site where the background 
electron density is pi. The embedding energy is given in 
the form 

^.(p-.) = ^^'''°'"'^"^^'^ ^''"'-°' (A3) 

\-A,E^^p, ifp, <0, 

where the sublimation energy _E? and parameter Ai de- 
pend on the element type of atom i. The background 
electron density pi is given by 



.(0) 



VG(r,;), 



(A4) 
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where 



r. = ^t 




(A5) 
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and 



G(r) = 



a/T 



-VW 



if r> -1, 
if r < -1. 



(A6) 



The zeroth and higher order densities, p\ \ p\ t Pi , and 

p] are given in Eqs. (A9|. The composition-dependent 
electron density scahng p^ is given by 



pO = p.oZ.oG (Tf) 



(A7) 



where pio is an element-dependent density scaling, Zio is 
the first nearest-neighbor coordination of the reference 
system, and FJ''^ is given by 



r 



rcf 



y2 Z-^ ''i "i ' 



(A8) 



fe=i 



Jk) ■ 



where s^ is the shape factor that depends on the ref- 
erence structure for atom i. Shape factors for various 



structures are specified in the work of Baskes' 
partial electron densities are given by 



m 



p. 



(0) 



Up. 



a(0) 



Tij) ^ij 



a 



T."' 



a(l) ^ija 



Oi 



j#« 



The 
(A9a) 

(A9b) 



1 

3 



a(2) '^ija'^ijfi 



P] ^ '^'' 



j¥=i 



E 



A"^) 



Vij ) ^t 



JT^i 



(A9c) 



(.?')'- E 



Q,/3,7 

3 



Ea(3) '^ija'^ij^'^ij'y 
Pj l3 



Oi. 



J/* 



E 






P. i>: 



J^i 



(A9d) 



where rij^ is the a component of the displacement vector 
from atom i to atom j. Sij is the screening function 



between atoms i and j and is defined in Eqs. (Al6l. The 
atomic electron densities are computed as 



a{k) I \ 

Pr [r,j) =p,oexp 



-#M^-1 



(AlO) 



where r,^ is the nearest-neighbor distance in the single- 

(k) 

element reference structure and (31 is element-dependent 
parameter. Finally, the average weighting factors are 
given by 



x(fe) 



(fc) MO) , 






(All) 



where ig ■ is an element-dependent parameter. 
The pair potential is given by 



4>ij inj) = 4>ij inj)S, 



J V ij ) ^i] 



(f>tj inj) 



Zji. 



-F, 



2i?^ iJ^^, 



'H] a(0) 



pr"(J^i) 



Z:. 



''-P^'^'in.) 



E^(r,,) = -E,,il + a*^(r,,))^<^^^'^^ 



a. 



*J 1 ^0 

2j 



(A12) 

(A13) 

(AM) 
(A15) 



where Eij, aij and r^ are element-dependent parameters 
and Zij depends upon the structure of the refe rence sys- 
tem. The background densities Pi{rij) in Eq. (A13) are 
the densities for the reference structure computed with 
interatomic spacing rij. 

The screening function Sij is designed so that Sij = 1 
if atoms i and j are unscreened and within the cutoff 
radius re, and 5^- = if they are completely screened or 
outside the cutoff radius. It varies smoothly between 
and 1 for partial screening. The total screening function 
is the product of a radial cutoff function and three body 
terms involving all other atoms in the system: 






•Jij — 'Jijjc 






^ikj — Jc 



Cikj — 1 + 2 



ij ik ~^ ij jk ij 



{r^k 



1 X>1 

fc{x)^\ [l-(l-a;)4)]' 0<x<l 
x<Q 



(A16a) 
(Al6b) 

(A16c) 
(A16d) 

(Al6e) 



Note that Cmin and C,„ax can be defined separately for 
each i-j-k triplet, based on their element types. The 
parameter Ar controls the distance over which the radial 
cutoff is smoothed from 1 to near r — re. 



Appendix B: Equilibrium lattice parameter and bulk 
modulus 

M E AM postulates the Rose universal equation of 
statePSi 



Ej,{a*) = -EAl + a* +5 



aa 



a + a* 



(Bl) 



for the reference structure of each single element and 
for each element pair. The a*, scaled distance from the 
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equilibrium nearest neighbor position tq, is 

a* — air/rQ — 1). 



(B2) 



Two 5 parameters may be specified for each element/pair: 
5r for negative, and 5a for positive a* . Then 



5r for a* < 
5a for a* > 0. 



(B3) 



The MEAM potential parameter a is related to the equi- 
librium atomic volume JIq, the bulk modulus i?o, and the 
cohesive energy of the reference structure E^ as follows 



QSoSlo 



Er. 



(B4) 



Appendix C: Trigonal and tetragonal shear modulus 

For small deformations of a cubic crystal, the change 
of energy density due to straining is 



^Ev = -Cu (e? + tl + el) + C12 (£162 + £263 + £3£i) 



+ -Cii{4 + tl + el)+0{4), 



(CI) 



where €i are strains in modified Voigt notation. 

The trigonal shear modulus C44 was determined from 
rhombohedral deformation given by £1 = £2 = £3 = and 
£4 = £5 = £6 = (5 in|Cl[ leading to 



(C2) 



A£;y(5) = -C44<5'+0((53). 



The DFT equilibrium energies and bulk moduli were 
obtained by fitting ener gy-v olume dependence to Mur- 
naghan equation of state^^Zl 



E{V)^E{V^) 

BqV 

^ B',{B',-l) 



B' I 



Vo 
V 



(B5) 



-1 



The tetragonal shear modulus (Cu — Ci2)/2 was deter- 

1 

1+5 



mined from the deformation given by £1 = (5, £2 = tts — 1 



CI leading to 

AEv{5) = {Cii-Ci2)5'' + 0{5^). (C3) 



The |C2| and |C3| were used to estimate tetragonal and 
trigonal shear moduli. 
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